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Abstract 

In toroidal geometry, and prior to the establishment of a fully developed turbulent state, the 
so-called topological instability of the pressure-gradient-driven turbulence is observed. In this 
intermediate state, a narrow spectral band of modes dominates the dynamics, giving rise to the 
formation of iso-surfaces of electric potential with a complicated topology. Since E x B advection 
of tracer particles takes place along these iso-surfaces, their topological complexity affects the 
characteristic features of radial and poloidal transport dramatically. In particular, they both 
become strongly non-diffusive and non-Gaussian. Since radial transport determines the system 
confinement properties and poloidal transport controls the equilibration dynamics (on any magnetic 
surface), the development of non-diffusive models in both directions is thus of physical interest. 
In previous work, a fractional model to describe radial transport was constructed by the authors. 
In this contribution, recent results on periodic fractional models are exploited for the construction 
of an effective model of poloidal transport. Numerical computations using a three-dimensional 
reduced magnetohydrodynamic set of equations are compared with analytical solutions of the 
fractional periodic model. It is shown that the aforementioned analytical solutions accurately 
describe poloidal transport, which turns out to be super diffusive with index a = 1. 
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I. INTRODUCTION 

Numerical calculations of resistive pressure-gradient-driven turbulence [l| in toroidal ge- 
ometry show the existence of an unstable regime below the threshold for fully developed 
turbulence. In this regime, the toroidal mode spectrum in steady state is dominated by a 
single toroidal mode N. This dominant mode may fluctuate intermittently among a narrow 
range of possible values. For instance, with the parameters used in Ref. [l], 24 < N < 26. 
The transition from a stable plasma to this ballooning-mode-dominated regime has the char- 
acteristic properties of a topological instability 0, [J. After the transition, the iso-surfaces 
of electrostatic potential induced by the ballooning modes have a complicated topological 
structure, which is a direct consequence of the (inner/outer) asymmetry in magnetic field 
strength inherent to any magnetic toroidal geometry. Consequently, particles advected by 
the E x B flows associated to these complex iso-surfaces cease to behave diffusively. At this 
point, it is worth stressing that this situation is quite different from cases of near-critical 
turbulence discussed in the literature 4j, in which the non-diffusive nature of transport is 
related instead to the existence of spatio-temporal correlations between fluctuations and the 
background profile gradients from which they feed. 

The structure of the potential iso-surfaces can be visualized as we move in the toroidal 
direction around the torus, following the magnetic field lines. At the singular magnetic sur- 
faces, (potential) vortices emerge with a structure which is consistent with the local twist of 
the magnetic field lines. As the outermost part of torus (the low-field side) is approached, 
filamentary vortices from different singular surfaces may merge and form extended radial 
streamers. Radial transport takes place predominantly within these streamers, since par- 
ticles can freely travel along them in the radial direction. The nature of radial transport 
within these streamers characterizes the confinement properties of the system. As we con- 
tinue to follow the lines toroidally, the vortices move back into the high-field side and the 
streamers break up. Each particle will then remain trapped within one of the filaments which 
emerge from this process. Poloidal (and toroidal) transport then follows, as the population 
of particles spreads out poloidally due to the free (ballistic) motion of each particle along 
the radially localized filament which contains it. Its nature thus determines how efficiently 
gradients are equilibrated on any magnetic surface. Eventually, some of the filaments will 
reach again the low-field side and merge to produce a new radial streamer, enabling again 
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radial transport until the structure goes back in to the high-field side, brakes again up into 
new filaments, and ballistic poloidal spreading ensues again. And so on. 

In Ref. jsj the radial transport of tracer particles in the unstable regime of this system 
was studied and found to be of a fractional, self-similar nature. This means that their trans- 
port can be described in terms of continuous-time random walks or fractional differential 
equations with appropriate exponents. It was also proven that the dynamics of tracers is 
pseudochaotic by showing that the dispersion of trajectories of (initially close) particles is 
polynomial (which in particular implies that the Lyapunov exponent is zero). In addition, an 
analogy between the topological structure of the flows and a billiard model (a paradigmatic 
example for pseudochaos) was discussed at length. 

In this paper we focus instead on the construction of an effective transport model that 
captures the main features of poloidal transport. This task is much more than a simple 
academic exercise or a straightforward extension of the radial transport case. The reason 
is that the difference with the radial case is not only the different dominant physics, but 
the existence of a periodic boundary. In the radial case, an absorbing boundary exists 
at the plasma edge, which can be dealt with within the standard framework of fractional 
differential operators. It has not been until very recently that the mathematical framework 
necessary to deal with a periodic boundary has been worked out by deriving the fluid limit 
equations of a Continuous-Time Random Walk formulated on a circle k3|. In the numerical 
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calculations of Ref. [11], the safety factor q was taken to be between 1 (at the plasma center) 
and 2 (at the plasma edge). Therefore, for a dominant toroidal mode N, there are N + 1 
possible filaments forming and the poloidal spreading of the tracers will take place in up to 
iV finite-size steps in q. As N increases, the poloidal spreading of particle tracers should 
approach the results for a 'fluid' transport model. In this paper, we will show that the 
new framework can deal successfully with the problem of poloidal transport in this system. 
It will also serve as illustration for future applications regarding transport along periodic 
directions in any magnetic configuration. 

The rest of this paper is organized as follows. Section [III gives a survey of the main results 
of on Continuous Time Random Walks and the fractional diffusion equation formulated 
on the circle. In Section imi we present the reduced magnetohydrodynamic model used in the 
study of the resistive pressure-gradient-driven turbulence and compare the numerical and 
analytical results for the transport of tracer particles. The conclusions are given in Section 
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IIVI The appendix contains some basic definitions on stable Levy distributions. 



II. FRACTIONAL DIFFUSION EQUATION ON A CIRCLE 

Continuous Time Random Walks (CTRWs) 6|, |7( are models describing the microscopic 
transport of particles in a probabilistic way. In this paper we are interested in the inter- 
pretation of poloidal transport in fusion plasmas as a CTRW defined on a circle. A general 
treatment of CTRWs on the circle has recently appeared In this section we collect the 
results of which are relevant for the present work. 

We denote by n(8,t) the density of particles (tracer particles for the application of 
the formalism relevant to this paper) normalized to the total number of particles, i.e. 
Jq W n(6,t)d8 = 1, Vt. The function n(9,t) must be periodic in 9, n{9 + 2ir,t) = n(9,t). 
A separable, Markovian, homogeneous time-translational invariant CTRW is defined by a 
mean waiting time, r, and a step-size pdf, p(A), giving the probability that a particle per- 
forms a jump from x to x+ A. The conservation of probability requires that J"_^p(A)dA = 1. 
Since A runs over the interval (— oo, oo), particles can wind around the circle an arbitrary 
number of times in each jump. We are interested in studying the case in which p(A) is a 
Levy stable distribution (see Appendix |Aj) . When the index of stability a = 2, the pdf p(A) 
is Gaussian, whereas for a < 2 it has algebraic tails. One would expect that (at least) in 
the latter case the effect of the non-trivial topology of the circle be very relevant and the 
CTRW on the circle should exhibit significant differences with respect to a CTRW with the 
same step-size pdf formulated on the real line. 

The dynamics described by the CTRW defined above is equivalent to the following Gen- 
eralized Master Equation (GME) [||: 

d t n(9, t) = - 1^ p(9 - 9')n(9', t')d9' - (1) 
r Jo T 

with 

oo 

p(0) = p{9 + 2nm). (2) 



m=— oo 



The sum in 



Notably, p is obtained from p by means of a ballooning transform 
explicitly accounts for the aforementioned fact that, given 9, 9' e [0, 2ir), particles can arrive 
at 9 from 9' through jumps of length \9 — 9' + 27rm|, m G Z. 
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Consider the case of a step-size pdf given by a stable Levy distribution (see Appendix 
with (3 — 0, whose characteristic function is: 



p(k) = exp(-a a \k\ a + ifik). (3) 

Then, the fluid limit of the GME (TTj) is: 

dtn = ^—nzW + 2 ^> + -d e n, (4) 

2r cos(7ra/z) r 

where o^ Q and 2n T> a are the Riemann-Liouville operators on the circle derived in Ref. 
The solution of Eq. (jl]) with initial condition n(9, 0) = X]-m=-oo^(^ — 27rm) (i.e. the 
propagator) is given by 



n(9,t) = ^- V e ( - CTQHQ+ ^ m) * /r e- ime , « G (0,2]. (5) 

27T Z — ' 

m=—oo 

In the limit t — > oo, n(9,t) — > 1/27T, as required by the conservation of the number of 
particles. 

When a = 1, which as we will see is the relevant case for the present work, the infinite 
sum on the right-hand side of Eq. (JHJ) can be computed and a closed expression for the 
propagator can be obtained: 

, n . 1 sinh(crt/r) , . 

Ltx cosn(o"r/rJ — cos(/it/r — 9) 

III. DYNAMICAL MODEL AND NUMERICAL CALCULATIONS OF TRACER 
PARTICLE TRANSPORT 



In toroidal geometry, the underlying instability of the resistive pressure-gradient-driven 
turbulence is the so-called resistive ballooning mode To calculate the dynamical prop- 
erties of these instabilities, we use the reduced set of Magnetohydrodynamics equations 
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ll| in toroidal geometry, and because of the low (3 values, we also use the electrostatic 
approximation. The model is then reduced to two equations: the perpendicular momentum 
balance equation and the equation of state. The former can be written in terms of the 
toroidal component of the vorticity U. In dimensionless form it reads: 

f = - S *B • V ( B • V*) + • Vp + ,VlU. (7, 



As for the equation of state: 



dp 



dt 



D\\— B Vf — B Vpj +D ± V 2 ±P . 




(8) 



Here, d/dt = d/dt + Vj_ • V is the convective derivative, and 



Vl 



— V$ x b 



(9) 



where $ is the electrostatic potential. There is a simple relation between the toroidal 
component of the vorticity and the stream function: 



Therefore, we solve Eqs. (J7J) and (jSJ) for p and <3>, taking into account Eqs. © and (|T0|) 
to relate U and $. In Eq. ([7J), 77 is the plasma resistivity, b = J5/B is a unit vector in the 
direction of the magnetic field, and re = b • Vb is the magnetic field line curvature. The 
magnetic field is expressed as B = FV( + V£ x Vip, where F = RB^ is the toroidal flux 
function, which is a very slowly varying function of the radial coordinate, r; ip is the poloidal 
flux, which is not evolved in time, and ( is the toroidal angle. Apart from the dissipation 
terms, we have two dimensionless parameters in these equations, j3 = p(0)/(B^/2p: ) and 
S = r R /Th p , the Lundquist number. Here, r R is the resistive time at the magnetic axis, 
t r — A i o a2 / ? ?(0), and r hp is the poloidal Alfven time, 77^ = R 0y /JIoff^n^/ B^, where rrii is 
the ion mass, and a and Ro are the minor and major radius, respectively. In the above 
dynamical equations lengths are normalized to the minor radius a, and time to the resistive 
time tr. 

To study the particle transport properties induced by these flow structures, we use pseudo- 
particles as tracers. These tracers are solutions of the equation of motion: 



Here, the velocity is the flow velocity given by Eq. @ in terms of the stream function, Vq is 
an arbitrary velocity along the field line, and b = £ — 6/q. Since this model is electrostatic, 
all information on turbulence evolution comes through the electrostatic potential $. In the 
present calculations we will use a constant Vq as initial condition for the tracer particles and 
keep the velocity field frozen in time because we are only looking for the effect of the flow 
structure on the transport. 



U = -C ■ V x Vi 



(10) 



dr 

dt 



Vi (r,t) + V b 



(11) 
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The first results of particle tracer transport have been obtained by launching 50000 
particles all with the same fixed initial velocity Vq = 4007T. At t — particles were located in 
a small region around r = 0.7, 8 = 0, and ( = 0. The evolution is initially very asymmetric 
because of the preferential direction induced by the drift in the motion of the particles along 
the eddies, which are aligned with the field lines. However, as particles go several times 
around the torus, n(8, t) becomes increasingly symmetric. For this reason, we first remove 
the drift of the distribution of tracers and then we compare the evolution with the symmetric 
solution of the fractional diffusion equation on the circle, Eq. 

The Lyapunov exponent in the motion of the particle tracers is zero and the trajectories 
separate from each other approximately linearly in time. Therefore, we expect a to be 
close to 1. Estimates based on simplified models give a G [1.05, 1.11] [j], [3]. To test this 
assumption and to determine the value of a associated with the poloidal transport, we first 
look at the time evolution of the width of the distribution, i.e. the square root of the second 
moment of the particle distribution, W(t). As shown in Fig. HJ W(t) increases linearly with 
time until it saturates at a constant value, W(oo) = vt/a/3, consistent with the flat tracer 
distribution. This behavior was expected from the analytical calculation. A fit of the data 
gives a = 0.984 ± 0.027, so that we can set a = 1. For the sake of completeness, we have 
plotted in Fig. [5] the time evolution of the first and third moments of the numerical pdf of 
tracers (after removing the drift), which are essentially zero at any time. 

According to Eq. ([6]), the expected time evolution of the peak of the particle distribution 
is given by 

11 + e~ at/T 

n(0,t) = — —r . (12) 

By fitting this expression to the numerical data, we can determine the time decay constant 
a/r. The result of the fit is shown in Fig. [31 and the value of the constant is a/r = 
0.268 ± 0.013. Now, we can use the analytical expression for n(6,t) given by Eq. ([6]) and 
compare this analytical prediction with the numerical data. We show a few examples of 
this comparison in Fig. HI The analytical solution seems to describe the numerical results 
relatively well. Recall that Eq. ( 1T21) is derived from a CTRW with symmetric step-size pdf, 
hence the numerical distribution of tracers should have vanishing odd moments, as is indeed 
the case (see Fig. [2]). 

We have also calculated the time evolution of particle tracers with random initial veloc- 
ities in the toroidal direction. We have used 50000 particles and the results of a typical 



calculation showing the evolution of n(6, t) are plotted in Fig. [5J In Fig. 0, one notes that 
the structures are gone. This is not surprising and it is one of the expected consequences 
of the random velocity initialization of the tracers. What may be surprising is the change 
of the functional form of the distribution. These distributions do not look at all like the 
analytical distributions given by Eq. ((6]). 

The explanation is simple. With a single initial velocity for all particles, at any given time 
they were all located at the same toroidal angle. Now, with the random initial velocities, at 
any given time the particles are distributed over a range of toroidal angles and what we have 
measured is a poloidal distribution averaged over all the toroidal angles. Furthermore, since 
there is a mean drift associated with the twist of the magnetic field lines, in each toroidal 
plane £, the poloidal distribution has its peak at a different value of 9. In each toroidal 
plane (, the analytical distribution is given by Eq. (JU]) with values a given value of ( in the 
range [— V max t, V max t\. Here V mscK = 4007T is the maximum velocity of the tracers. Therefore, 
to reproduce the measured poloidal distribution, we have to average the fixed ( distribution 
over the toroidal angle, that is 

1 - e- 2 ^ 1 f d( 

Wl ' } 2tt 2V max t J l-2e-<^cos(0-O + e-^' 1 J 

Vmax^ 

where u is the average pitch of the field line at the radial position of the initial tracers. After 
some algebra, one obtains 



ln)(9,t) = - 



1 f ri + e ^'/ r / 

< arctan — tan 



+ uV max t 



— arctan 



1 + ^ (Q- uV mKK t 
-^tan 

1 - e ~ at l T 



(14) 



This distribution should correspond to the one obtained in the numerical calculations. 
Using the value of a obtained from the previous section and the following estimate for 
MV max , u = 2/3 and V max = 4007T as used in the numerical calculations, we have plotted 
the corresponding pdfs in Fig. [Ul The agreement between Figs. and E] is very good. For 
the randomized velocity case, once we have averaged over the initial velocities, there is no 
difference between anomalous diffusion with a = 1 and pure ballistic motion of the particles 
in the absence of the flow structures. However, the distinction is clear when we look at 
monoenergetic particles, as we saw at the beginning of this section. 
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IV. CONCLUSIONS 



In this paper, we have used numerical simulations of resistive ballooning mode turbulence 
at low beta (sufficiently so to produce a topological instability) to illustrate the adequacy 
of recent periodic formulations of fractional transport equations Q to describe non-diffusive 
turbulent transport in the poloidal direction. The formalism is quite general and could be 
easily applied to many other instances of transport problems along a periodic direction, of 
which parallel equilibration dynamics in toroidal magnetic devices is just one example. 

In the case examined here, non-diffusive transport ensues due to the complex topology of 
the potential iso-surfaces present in the system, a direct consequence of the narrow spectral 
band of dominant modes. Both radial and poloidal transport exhibit non-diffusive features 
due to this complex topology, but their physical origin is rather different, as explained in 
the text. The radial anomalous diffusion in this system was already considered in Ref. 13]. 
In the current paper, we have focused instead on the poloidal transport, which turns out to 
be an interesting example of a = 1 anomalous diffusion. The origin of this value of a is in 
the fact that, away from the streamers, particles simply spread out ballistically along the 
filaments while on the high-field side. Once they reach the low-field side again and enter 
a new streamer structure, their direction of motion can be modified as the new streamer 
breaks up to form a new set of filaments towards the high-field side and particles get trapped 
in them. It is this combination of ballistic motion plus scattering of velocity directions which 
is ultimately responsible of a Cauchy type (i.e., a = 1) diffusive process. 

By looking in detail at the poloidal distribution of monoenergetic particle tracers, we 
have also shown that the new mathematical framework can capture the features of the 
poloidal transport process. There is very good agreement between numerical results and the 
analytical calculation describing this phenomenon (sf] . If tracers with random velocities are 
used, once we average over initial conditions, it is not possible to distinguish the anomalous 
diffusion with a = 1 from the normal ballistic motion of the particles in the absence of the 
flow structures. 
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Appendix A: LEVY SKEW ALPHA-STABLE DISTRIBUTIONS 

The family of Levy skew alpha-stable distributions (or simply stable distributions, or Levy 
distributions) is parameterized by four real numbers a G (0,2], f3 G [—1,1], a > 0, and 
Their characteristic function (i.e. their Fourier transform) is given by [jjj]: 

exp (-a a \k\ a [l - i/3sign(fc) tan(^)] + ifik) a^l, 

(Al) 

exp (— <r|A;| [l + z/5^sign(/c)ln|/c|j + ifik) a — 1. 



According to the Generalized Central Limit Theorem [14]], stable distributions are the 
only possible distributions with a domain of attraction. The index a is related to the 
asymptotic behaviour of S(a, (3,a, fi)(x) at large x: 

[ C <* (^) ^M" 1 "" x -> -oo, 
S(a,P,<T,fi)(x)= I (A2) 

{ C a (^) X^OO, 

for a G (0,2). For a = 2, S(2, j3, a, /i) is a Gaussian distribution. 
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Figure 1: Standard deviation of the pdf of particle tracers versus time. The short-time behaviour 
implies that the index a is equal to 1. 
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Figure 3: Numerical results for n(0, t) obtained from Eq. (jlip are fitted to the analytical expression 
Eq. CE}, yielding a/r = 0.268 ± 0.013. 
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Figure 4: Comparison of the time evolution of the pdf of tracer particles with fixed initial toroidal 
velocity Vq = 4007T, and the analytical expression, Eq. 0. 
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Figure 5: Numerical computation of the time evolution of the pdf of tracer particles with random 
initial toroidal velocities. 
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Figure 6: Time evolution of n(6, t) given by Eq. (|14p for the values of the parameters corresponding 
to the case plotted in Fig. 
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